Phase transitions in systems of magnetic dipoles on a square lattice with quenched 
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We study by Monte Carlo simulations the effect of quenched orientational disorder in systems of 
interacting classical dipoles on a square lattice. Each dipole can lie along any of two perpendicular 
axes that form an angle ?/) with the principal axes of the lattice. We choose tp at random and 
without bias from the interval [—A, A] for each site of the lattice. For < A < 7r/4 we find a 
thermally driven second order transition between a paramagnetic and a dipolar antiferromagnetic 
order phase and critical exponents that change continously with A. Near the case of maximum 
disorder A « 7r/4 we still find a second order transition at a finite temperature Tc but our results 
point to weak instead of strong long-ranged dipolar order for temperatures below Tc. 

PACS numbers: 75.10.-b, 75.40.Cx, 75.40.Mg 
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I. INTRODUCTION 

Systems of interacting dipoles (SIDs) are atracting a 
renewed interest. This is in part due to recent advances in 
nanoscience^ which render available realizations of mag- 
netic nanoparticle assemblies^. These systems show a 
rich collective behaviour in which the dipole-dipole in- 
teractions play an essential role. This is because dipole- 
dipole interaction strength grows linearly with nanoparti- 
cle volume. Collective effects that are controlled by dipo- 
lar interactions can therefore become important at rea- 
sonably high temperatures. Spatial variations of the di- 
rection of magnetic dipolar fields lead to frustration and 
make SIDs very sensitive to their spatial arrangements>3.. 
In crystalline arrangements for example, SIDs exhibit 
different long ranged ferro or antiferromagnetic dipolar 
magnetic order that depends crucially on lattice geom- 
etry. The magnetic ordering depends also crucially in 
anisotropy. On the one hand, single site anisotropy is 
always present in magnetic crystals. On the other hand, 
dipole-dipole interactions by themselves create effective 
anisotropics in 2D systems. In square lattices, for ex- 
ample, dipole-dipole interactions pushes spins to lie on 
the plane of the lattice, and termal excitations tends to 
align them along the two principal axes of the lattice^. 
The resulting magnetic order states that ensue from the 
competition of dipolar and anisotropic energies could be 
as exotic as "spin ice" found in diamond type crystals^. 
SIDs in disordered spatial arrangements are also of in- 
terest. Some non equilibrium spin glass behaviour (like 
time dependent susceptibilities and aging) has been ob- 
served in experiments with systems of randomly placed 
dipoles as frozen ferrofluids or strongly diluted magnetic 
crystals^. Very recent computer simulations indicate that 
systems of Ising dipoles with random anisotroy axes have 
an equilibrium spin glass phase at low temperatures^. 

The aim of this paper is to study by numerical sim- 
ulations the effect of quenched directional disorder in a 
system of dipoles placed in cristalline array. We are in- 



terested in how the magnetic dipolar order and the char- 
acter of the transition between the paramagnetic and the 
low temperature ordered phase vary as we increase the 
amount of disorder. We also investigate whether there is 
a well defined threshold of disorder beyond which dipolar 
order dissapears. 



II. MODEL AND CALCULATION 

We first define the system model we use. Let Si be a 
classical two xy-component unit spin at lattice site i of 
a square lattice. These spins interact as pure magnetic 
dipoles with Hamiltonian 
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FIG. 1: (a) Order parameter vs T for systems with A = 
0.4. X, A, o, □ and o stand for L = 8,12,26,22 and 32 
respectively. Lines are only guides to the eyes. In the inset, 
specific heat vs. T for the same A and systems sizes, (b) and 
inset therein, same as in (a) but for A = 7r/4. 
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Yij is the displacement from site i to site j, a is the square 
lattice parameter. In the following all energies and tem- 
peratures are given in terms of and ed/ks respectively. 
When there is no quenched disorder, we consider a strong 
quadrupolar anisotropy that forces spins to lie along any 
of the two principal axis of the square lattice. This 
dipolar four-state clock model has been studied recently, 
and found to have a second order transition between a 
paramagnetic and an antiferromagnetic dipolar phase at 
Tc = 1.106 with critical exponents a/v = 0.82(5) and 
P/v — 0.03(2)5. In the ordered phase spins point up 
along lines with alternate sign from one line to the adja- 
cent one. In order to characterize these antiferromagnetic 
order, we define a staggered magnetization as 



(3) 



and use (j) = ('t^x + ti^)^/^ as order parameter. 

This system has been found to be very sensitive to 
the direction of the quadrupolar anisotropy axes. For 
example, if one turn these axes an angle ip = tt/A, critical 
exponents change to a/iy = 0.44(4) and P/v = 0.00(5)^. 
Thus, we introduce disorder in the model described above 
by tilting these two orthogonal anisotropy axes by an 
angle ip with respect to the crystalline axes. For each spin 
this angle V' is chosen at random and without bias from 
the interval [—A, A]. Note that A — 7r/4 correspond 
to a completely random distribution of the orientation 
of these pairs of axes. Note also that even in this case 
anisotropics never force any pair of spins to form angles 
greater that 7r/4. Therefore, this is less disruptive than 
random uniaxial anisotropy^. 

We study systems of L x L spins and use periodic 
boundary conditions. We let a spin on site i interact 
through dipolar fields with all spins within an L x L 
square centered on the i- site*. Our simulations follow 
the standard Metropolis Monte Carlo (MC) algorithmjiS. 
We use a single-spin flip dynamics, in which all dipolar 
flelds are updated throughout the system every time a 
spin flip is accepted, before another spin is chosen in or- 
der to repeat the process. In our MC simulations, we 
start for a given realization of disorder from T — 2.5 
(well in the paramagnetic phase) and lower temperature 
in AT = 0.05 steps. At each value of T we averaged 
over 5 x 10^ MC sweeps having flrst discarded 10^ MC 
sweeps to let the system equilibrate. Finally, our results 
are averaged over Nr different independent realizations 
of quenched disorder. In order to get realiable results 
we use at least Nr — 200 for A > 0.6 and Nr — 50 for 
A < 0.6. As a result, error bars in the figures shown 
in this paper are always smaller than symbol sizes used 
therein. Our results follow from simulations for system 
sizes L = 8, 12, 16, 22 and 32. This would be insuficient to 



obtain accurate critical exponent values but is adequate 
for our purposes. 

We study the modulus of the staggered magnetiza- 
tion (j) defined above and obtain the susceptibility from 
their fluctuations as x = ScfP' / {NT). We calculate the 
speciflc heat from energy fluctuations via the relation 
C = SE"^ / [NT"^). Finally, we use the cumulant like quan- 
tity 



(4) 



As L ^ oo we expect that ui2 ^ in the paramagnetic 
phase, ui2 ^ 1 in the long ranged ordered phase and tend 
to some intermediate value at a critical point. Hence, 
curves of ui2 versus T for different system sizes should 
cross at critical points for large enough L. 



III. RESULTS 

Data obtained for the order parameter for A = 0.4 
are exhibited in a </> versus T plot in figure 1(a). Curves 
for different L cross at Tc = 0.95(1). Below this temper- 
ature <j) increases as L increases indicating the existence 
of an ordered phase. Plots for A = 0.2, 0.4, 0.65 and 
0.7 (not shown) are qualitalively similar, except in that 
curves for L > 16 tend to collapse below T^ for the case 
A = 0.7. The behavior is markedly different for A k. 7r/4 
(see Fig. 1(b)), where (j) appears to decrease as L increases 
even at low T at least for the sizes we have studied, rais- 
ing the question over the existence of strong long-range 
order. We return to this point in the discussion of Fig. 4. 

In the inset of Fig. 1(a) we plot the specific heat C 
vs T for A — 0.4. Clearly a singularity develops as L 
increases. We observe a weakening of this singularity 
as we increase A from to 0.7. For A k, 7r/4 we see 
that the singularity washes out completely (see inset of 
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FIG. 2: Log-log plot or the maximum value of the specific 
heat as a function of L for the indicated values of A. • stands 
for the model without disorder, but with quadrupolar easy 
axis forming angles of pi/4 with the x^y axis at every site of 
the lattice. 
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Fig. 1(b)). In Fig. 2 we show a log-log plot of Cmax 
versus T, where Cmax is the value of the specific heat 
at its maximum. We obtain a/u from the straight line 
slopes using Cmax oa L"/'^. This is so for A — 0.2, 0.4 and 
0.6 for which we obtain a/i^ — 0.83(4), 0.63(4) and and 
0.31(6) respectively. In contrast, data plotted in Fig. 2 
for A > 0.6 show not straight lines but curves that be- 
come flatter as L increases suggesting that a/v = 0. This 
dependence of a/v with A is in agreement with the Har- 
ris criterionii, in principle valid only for systems with 
short ranged interactions. 

We next discuss whether strong or weak long-range 
order exists at low temperatures for A = 7r/4. For this 
purpose, we compare log-log plots of vs 1/L for A ~ 0.6 
and A — 7r/4. These plots are shown in Figs. 3(a) and 
3(b) respectively and are qualitalively different. Fig. 3(a) 
is consistent with a phase transition from a magnetically 
disorder state above T = 0.95(1) to a strong long range 
order below them. We have obtained similar behaviour 
for other values of A in the (0, 0.7) range. On the other 
hand, in Fig. 3(b) appears to decrease algebraically 
with L over a wide range of temperatures for T < 0.68 
Dashed thick lines in the figure stand for linear regression 
of the data for T = 0.52 and 0.68 and their slope gives 
P/f = 0.06(2) and 0.12(4) respectively using </> oc L^^/". 

We also gathered information from the behaviour of 
susceptibility x- In Fig. 4(a), plots of x vs T for A = 7r/4 
show curves whose peak grows as L increases. From 
the data shown, we obtain Xmax oc L^-^^^^") where Xmax 
stands for the maximum value of x versus T for a given 
value of L. We obtain qualitatively similar results for 
A < 7r/4. Note also in Fig. 4(a) that the position of the 
maximum of x, Tm, changes with L. In the inset of Fig. 5 
we plot Tm versus 1/L for different values of A. Direct 
extrapolation of these data to 1/L —> gives an esti- 
mation of the transition temperature Tc. The resulting 
values of Tc vs A are plotted in Fig. 5. 
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FIG. 4: (a) Susceptibility versus T for A = 7r/4. x. A, o, 
□ and o stand for L = 8, 12, 26, 22 and 32 respectively Lines 
are only guides to the eyes, (b) Cumulant versus T for 
A = 7r/4. Symbols stand for the same system sizes as in 
(a). In the inset of Fig. (b), u\2 versus l/L for different 
temperatures and A = n/A:. o, •, ■, o, ▲ and x stand for 
T = 0.80, 0.76, 0.72, 0.68, 0.64 and 0.52 respectively 
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A more accurate determination of Tc and some ad- 
ditional information about the nature of the transition 
was obtained from ui2. In Fig. 4(a) we plot ui2 vs T 
for A = 7r/4. Curves for different values of L exhibit a 
clear crossing at = 0.69(1) and ui2 = 0.82(1). We 
have obtained similar plots for various values of A that 
enabled us to estimate Tj, by the value of T where the 
curves cross. The resulting values of Tc are plotted vs A 
in Fig. 5, which gives the global phase diagram for our 
model. In the inset of Fig. 4(b) we plot ui2 versus 1/L 
for A — 7r/2. For the system sizes we have considered, 
curves exhibit marked finite size effects even well below 
T = 0.69. Results for T = 0.52 and 0.64 shown in this in- 
set seem to indicate that ui2 ^ 1 as 1/L ^ or at least 
that 1<<0^ >/<(/) >^< 1.02. This is consistent with 
the behaviour— found for ui2 for 2D XY — like models 
that exhibit a Kosterlitz-Thouless transition with weak 
long-range order below a transition temperature. 



In sum, we have reported evidence from MC simula- 
tions that disorder in the orientation of the quadrupo- 
lar anisotropy axes on systems of interacting dipoles 
in square lattice is relevant, in the sense that modifies 
the critical behavior of the thermal transition between a 
paramagnetic and the antiferromagnetic dipolar phase. 



FIG. 3; (a) Log-log plot of the order parameter versus 1/L for 
A = 0.4 and the temperatures indicated in the figure. Error 
bars are smaller than the size of the symbols, (b) Same as 
in (a) but for A = 7r/4. Dashed thick lines in (b) stand for 
linear regression of the data. Thin solid lines are guides to 
the eyes. 
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FIG. 5: Transition temperature Tc versus A. Q stand for 
values of Tc obtained from extrapolation of temperatures at 
which susceptibility attains its maximum for different values 
of L (see text), x stand for Tc obtained as the temperature 
for which u\2 is independent of L. In the inset, temperature 
of maximum susceptibility versus 1/L. A, o, □ and o stand 
for A = 0.4, 0.6, 0.7 and 7r/2 respectively. 
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